
/*
*******************************************************************************
* ansi c source code
* file name:
*	Random.c
* abstract:
*	to provide a series of functions, which can offer random numbers
* author:Fangming He
*******************************************************************************
*/


/*
*******************************************************************************
*                               include files
*******************************************************************************
*/
#include <stdio.h>
#include <string.h>
#include <math.h>
#include "globaldef.h"


/*
*******************************************************************************
*			           constants and define declarations
*******************************************************************************
*/
/* for function RunRand1 */
#define IA1		16807
#define IM1		2147483647
#define AM1		(1.0/IM1)
#define IQ1		127773
#define IR1		2836
#define NTAB1	32
#define NDIV1	(1 + (IM1 - 1)/NTAB1)
#define EPS1	(1.2e-7)
#define RNMX1	(1.0 - EPS1)


/*
*******************************************************************************
* 		                    global object definition
*******************************************************************************
*/
/* local variables for function RunRand1() */
long	rand1_seed = 1;
long	rand1_iy = 0;
long	rand1_iv[NTAB1];


/* local variables for function RunRand2() */
long	rand2_seed = 1;
long	rand2_x;
long	rand2_y;
long	rand2_z;
int		rand2_cnt = 0;
/*
*******************************************************************************
*                            function definition
*******************************************************************************
*/
/*
*******************************************************************************
* description:
*	initialize the local variables for function RunRand1()
* input: 
*	seed: the seed for generating random figures
* output:
*	none
* function reference:
* author:
*	guo liang 
*******************************************************************************
*/
void SetRand1(long seed)
{
	rand1_seed = seed;
	rand1_iy = 0;
	memset(rand1_iv, 0, NTAB1 * sizeof(long));
}


/*
*******************************************************************************
* description:
*	to produce uniform iid sequence within [0.0 1.0)
* input: 
*	none
* output:
*	return: a uniformly distributed number within [0.0 1.0)
* function reference:
* author:
*	guo liang 
*******************************************************************************
*/

double RunRand1(void)
{
	int		j;
	long	k;
	double	temp;

	if ( (rand1_seed <= 0) || (!rand1_iy) )
	{
		if ( -rand1_seed < 1 ) rand1_seed = 1;
		else rand1_seed = -rand1_seed;

		for (j = NTAB1+7; j >= 0; j--)
		{
			k = rand1_seed / IQ1;
			rand1_seed = IA1 * (rand1_seed - k*IQ1) - IR1*k;

			if ( rand1_seed < 0 ) rand1_seed += IM1;
			if ( j < NTAB1 ) rand1_iv[j] = rand1_seed;
		}

		rand1_iy = rand1_iv[0];
	}

	k = rand1_seed / IQ1;
	rand1_seed = IA1 * (rand1_seed - k*IQ1) - IR1*k;
	if ( rand1_seed < 0 ) rand1_seed += IM1;
	j = rand1_iy / NDIV1;
	rand1_iy = rand1_iv[j];
	rand1_iv[j] = rand1_seed;
	temp = AM1 * rand1_iy;
	if ( temp > RNMX1 ) return RNMX1;
	else return temp;
}
/*
*******************************************************************************
* description:
*	initialize the local variables for function RunRand4()
* input: 
*	seed: the seed for generating random figures
* output:
*	none
* function reference:
* author:
*	guo liang
*******************************************************************************
*/
void SetRand2(long seed)
{
	rand2_seed = seed;
	rand2_x = 0;
	rand2_y = 0;
	rand2_z = 0;
	rand2_cnt = 0;
}


/*
*******************************************************************************
* description:
*	to produce uniform iid sequence within [0.0 1.0)
* input: 
*	none
* output:
*	return: a uniformly distributed number within [0.0 1.0)
* function reference:
* author:
*	guo liang
*******************************************************************************
*/
double RunRand2(void)
{
	double	result;
	
	if ( rand2_cnt == 0 )
	{
		rand2_x = rand2_seed;
		rand2_y = rand2_seed;
		rand2_z = rand2_seed;
		rand2_cnt++ ;
	}

	rand2_x = (rand2_x * 249) % 61967;
	rand2_y = (rand2_y * 251) % 63443;
	rand2_z = (rand2_z * 252) % 63599;
	result = (rand2_x/61967.0 + rand2_y/63443.0 + rand2_z/63599.0) - (int)(rand2_x/61967.0 + rand2_y/63443.0 + rand2_z/63599.0);

	return result;
}
/*
*******************************************************************************
* description:
*	initialize the local variables for function RunRandGas1()
* input: 
*	seed: the seed for generating random figures
* output:
*	none
* function reference:
* author:
*	guo liang 
*******************************************************************************
*/
void SetRandGas(long seed)
{
	SetRand1(seed);
  
}

/*
*******************************************************************************
* description:
*	to produce complex Gassiuan distributed sequence 
* input: 
*	none
* output:
*	return: a complex Gassiuan distributed number with corresponding mean and variance
* function reference:
* author:
*	guo liang
*******************************************************************************
*/
struct Complex RunRandGas1(double mean,double variance)
{
	struct Complex result;
	double v1,v2;
	double coeff;
	double tmp;
	do
	{
		v1 = 2 * RunRand1() - 1;
		v2 = 2*  RunRand1() - 1;
		tmp = v1*v1 + v2*v2;
	}
	while( tmp>=1 || tmp == 0);
    coeff = sqrt(-2.0 * log(tmp) / tmp);
    result.real = v1*coeff*sqrt(variance) + mean;
	result.imag = v2*coeff*sqrt(variance) + mean;
	return result;
}
